clear all

set seed 2334239
local bsiter = 1000



cd "/Users/sandroambuehl/Dropbox/Demand for Information/DATA_ANALYSIS/"



set more off

use "data/data_all.dta", replace

drop if session == 1

merge m:1 subj using "data/alphas.dta"
drop _merge


save "data/pool.dta", replace


gen low00 = p00 < 0.47
gen low11 = p11 < 0.47
gen low = max(low11,low00)
sort subj
by subj: egen rogue2num = sum(low)
gen rogue2 = rogue2num >= 3


* get data and predict bayes posterior at hit


keep if tag_sub == 1
keep sub id alpha se_alpha session group rogue2 knows_bayes gender lotr1 lotr2 lotr3 lotr4 lotr5 lotr6 superstition1 superstition2 superstition3 religious_beliefs1 religious_beliefs2 religious_beliefs3 major crt_score race
gen id_new = _n

expandby 4, by(id_new)
sort id_new
gen info_struct = _n
replace info_struct = info_struct - 4*floor(info_struct/4)
replace info_struct = 4 if info_struct == 0
gen period = info_struct + 20

gen id_inf = 10*id_new + info_st
expandby 12, by(id_inf)
sort id_inf
gen round = _n
replace round = round - 12*floor(round / 12)
replace round = 12 if round == 0

qui merge m:1 subject period round  using "data/hitting.dta"
drop if _merge == 2
drop _merge

local t = 1

gen pred_logit_post = alpha * log(b_post / (1-b_post))

gen pred_hitting_round = round if 1 - (pred_logit_post < logit_outside_chance & pred_logit_post  > - logit_outside_chance)
sort subject period
by subject period: egen pred_first = min(pred_hitting_round)
replace pred_hitting_round = . if pred_hitting_round != pred_first
replace pred_hitting_round = 1 if pred_hitting_round != .

sort subj info_s
by subj info_s: egen pred_first_2 = min(pred_hitting_round)
replace pred_first = 13 if pred_first_2 == .
by subj: egen max_pred_first = max(pred_first_2)
replace max_pred_first = 13 if max_pred_first == .

gen pred_logit_b_post_when_hit = log(b_post / (1 - b_post)) if pred_hitting_round == 1
replace pred_logit_b_post_when_hit = abs(pred_logit_b_post_when_hit) //ensure that this is the posterior for the box that would have been chosen!

gen diff_pred_logit_b_post_when_hit = pred_logit_b_post_when_hit - logit_outside_chance

egen tag_sub_per = tag(subject period)

sort subj
by subj: egen max_first = max(first)


rename diff_logit_b_post_when_hit diff
rename diff_pred_logit_b_post_when_hit diff_sim

sort subject
by subj: egen mean_diff_pred = mean(diff_sim)
by subj: egen mean_diff = mean(diff)
gen logit_b_post = log(b_post / (1-b_post))
sort subj info_s
by subj info_s: egen max_b_post = max(abs(logit_b_post))

collapse rogue2 se_alpha mean_diff mean_diff_pred diff diff_sim max_first max_pred_first alpha session max_b_post logit_outside_chance knows_bayes gender lotr1 lotr2 lotr3 lotr4 lotr5 lotr6 superstition1 superstition2 superstition3 religious_beliefs1 religious_beliefs2 religious_beliefs3 major crt_score race, by(subject info_s) 
egen tag_sub = tag(subj)







****************************************
* tobit analysis ***********************
****************************************


* generate predicted delta^b for all

gen loglik = .
replace loglik = .84729779 if info_s == 1
replace loglik = 1.3862945 if info_s == 2
replace loglik = .84729779 if info_s == 3
replace loglik = .40546522 if info_s == 4	 

local steps = 10

forvalues i = 1 / `steps' {
	gen b_`i' = `i'*loglik
	gen a_`i' = abs(alpha) * b_`i'
	gen nohit_`i' = a_`i' < logit_outside_chance	// nohit is 1 if the beliefs haven't yet been sufficiently extreme
}

gen sum_1 = loglik * nohit_1 * 2	// 
forvalues i = 2 / `steps' {
	local iminus = `i'-1
	gen sum_`i' = loglik* nohit_`i' + sum_`iminus'
}
rename sum_`steps'  b_post_at_hit
gen diff2 = b_post_at_hit - logit_outside_chance

// replace diff_sim = diff2 //if diff_sim == .


gen afram = race == 1
gen asian = race == 2
gen cauc = race == 3

gen censor_limit = max_b_post - logit_outside_chance
gen censored = diff == .
replace diff = -10 if censored == 1


program define progMLE
	args lnf gamma0 gamma1 theta
	tempvar x y bound censored
	quietly{ 
		gen double `y' = $ML_y1
		gen double `x' = $ML_y2
		gen double `censored' = $ML_y3
		gen double `bound' = $ML_y4
	
		replace `lnf' = log( normal(`gamma0'+`gamma1'*`x'-`theta'*`bound')) if `censored' == 1
		replace `lnf' = log(`theta') + log(normalden(`theta'*`y' - `gamma0' - `gamma1'*`x')) if `censored' == 0
	}
end



program define progMLEcontrols
	args lnf gamma0 gamma1 theta beta1 beta2 beta3 beta4
	tempvar x y z1 z2 z3 z4 bound censored
	quietly{ 
		gen double `y' = $ML_y1
		gen double `x' = $ML_y2
		gen double `censored' = $ML_y3
		gen double `bound' = $ML_y4
		gen double `z1' = $ML_y5
		gen double `z2' = $ML_y6
		gen double `z3' = $ML_y7
		gen double `z4' = $ML_y8

		replace `lnf' = log( normal(`gamma0'+`gamma1'*`x' + `beta1'*`z1' + `beta2'*`z2' + `beta3'*`z3' + `beta4'*`z4' - `theta'*`bound')) if `censored' == 1
		replace `lnf' = log(`theta') + log(normalden(`theta'*`y' - `gamma0' - `gamma1'*`x' - (`beta1'*`z1' + `beta2'*`z2' + `beta3'*`z3' + `beta4'*`z4'))) if `censored' == 0
	}
end





postfile gradInfoBootstrap id b_1 b_cons r2 equation using "data/gradInfoBootstrap", replace

forvalues d = 1/`bsiter' {


preserve

bsample, cluster(subj)


****** Column 1

gen predictor = alpha

ml model lf progMLE (gamma0: diff predictor censored censor_limit =  ) (gamma1: ) (theta: ), technique(bhhh) init(0 1 1,copy) vce(cluster subject)

capture: ml maximize

local gamma0 = [gamma0]_b[_cons]
local gamma1 = [gamma1]_b[_cons]
local theta = [theta]_b[_cons]
local sigma = 1/`theta'
local beta0 = `sigma'*`gamma0'
local beta1 = `sigma'*`gamma1'

gen diff_star = `beta0' + `beta1'*predictor
su diff_star
local v = r(Var)
local MZ_R2 = `v'/(`v'+(`sigma'^2))
di `MZ_R2'

post gradInfoBootstrap (`d') (`beta1') (`beta0') (`MZ_R2') (1) 

drop diff_star predictor

****** Column 2

gen predictor = alpha

ml model lf progMLEcontrols (gamma0: diff predictor censored censor_limit afram asian cauc gender =  ) (gamma1: ) (theta: ) (beta1: ) (beta2: ) (beta3: ) (beta4: ), technique(bhhh) init(0 1 1 1 1 1 1,copy) vce(cluster subject)

capture: ml maximize

local gamma0 = [gamma0]_b[_cons]
local gamma1 = [gamma1]_b[_cons]
local theta = [theta]_b[_cons]
local sigma = 1/`theta'
local beta0 = `sigma'*`gamma0'
local beta1 = `sigma'*`gamma1'

gen diff_star = `beta0' + `beta1'*predictor
su diff_star
local v = r(Var)
local MZ_R2 = `v'/(`v'+(`sigma'^2))
di `MZ_R2'

post gradInfoBootstrap (`d') (`beta1') (`beta0') (`MZ_R2') (2) 

drop diff_star predictor


****** Column 3

gen predictor = diff2

ml model lf progMLE (gamma0: diff predictor censored censor_limit =  ) (gamma1: ) (theta: ), technique(bhhh) init(0 1 1,copy) vce(cluster subject)

capture: ml maximize

local gamma0 = [gamma0]_b[_cons]
local gamma1 = [gamma1]_b[_cons]
local theta = [theta]_b[_cons]
local sigma = 1/`theta'
local beta0 = `sigma'*`gamma0'
local beta1 = `sigma'*`gamma1'

gen diff_star = `beta0' + `beta1'*predictor
su diff_star
local v = r(Var)
local MZ_R2 = `v'/(`v'+(`sigma'^2))
di `MZ_R2'

post gradInfoBootstrap (`d') (`beta1') (`beta0') (`MZ_R2') (3) 

drop diff_star predictor

****** Column 4

gen predictor = diff2

ml model lf progMLEcontrols (gamma0: diff predictor censored censor_limit afram asian cauc gender =  ) (gamma1: ) (theta: ) (beta1: ) (beta2: ) (beta3: ) (beta4: ), technique(bhhh) init(0 1 1 1 1 1 1,copy) vce(cluster subject)

capture: ml maximize

local gamma0 = [gamma0]_b[_cons]
local gamma1 = [gamma1]_b[_cons]
local theta = [theta]_b[_cons]
local sigma = 1/`theta'
local beta0 = `sigma'*`gamma0'
local beta1 = `sigma'*`gamma1'

gen diff_star = `beta0' + `beta1'*predictor
su diff_star
local v = r(Var)
local MZ_R2 = `v'/(`v'+(`sigma'^2))
di `MZ_R2'

drop diff_star predictor

post gradInfoBootstrap (`d') (`beta1') (`beta0') (`MZ_R2') (4) 

restore

}

postclose gradInfoBootstrap

clear

use "data/gradInfoBootstrap", replace
noisily: su b_1 	if equation == 1
noisily: su b_cons 	if equation == 1
noisily: su b_1 	if equation == 2
noisily: su b_cons 	if equation == 2
noisily: su b_1 	if equation == 3
noisily: su b_cons 	if equation == 3
noisily: su b_1 	if equation == 4
noisily: su b_cons 	if equation == 4




